A massive helium star with a sufficiently strong 
magnetic field to form a magnetar 


Tomer Shenar!, Gregg Wade?, Pablo Marchant?, Stefano Bagnulo‘ 

Julia Bodensteiner??, Dominic M. Bowman, Avishai Gilkis?, Norbert Langer”, 
André Nicolas-Chené?, Lidia Oskinova!, Timothy Van Reeth?, Hugues Sana?, 
Nicole St-Louis!!, Alexandre Soares de Oliveira!?, Helge Todt!?, Silvia Toonen! 
lAnton Pannekoek Institute for Astronomy, University of Amsterdam, Amsterdam 1098 XH , the Netherlands 
? Department of Physics and Space Science, Royal Military College of Canada, Kingston K7K7B4, Canada 
3 Institute of Astronomy, Katholieke Universiteit (KU) Leuven, Leuven 3001, Belgium 
^ Armagh Observatory & Planetarium, College Hill, Armagh BT61 9DG, UK 
5 European Southern Observatory, Garching bei München 85748, Germany 
6 The School of Physics and Astronomy, Tel Aviv University, Tel Aviv 6997801, Israel 
7 Argelander-Institut für Astronomie, Universität Bonn, Bonn 53121, Germany 
8 Max-Planck-Institut für Radioastronomie, Bonn 53121, Germany 
9 National Science Foundation’s National Optical-Infrared Astronomy Research Laboratory (NSF's NOIRLab), Hawaiʻi 96720, USA 
10 Institut für Physik und Astronomie, Universitat Potsdam, Potsdam D-14476, Germany 
11 Département de physique, Université de Montréal, Complexe des sciences, Montréal H2V 0B3, Canada 
12 Institute of Research and Development, Universidade do Vale do Paraíba, 12244-000, Sáo José dos Campos, Brazil 


*To whom correspondence should be addressed; E-mail: T.Shenar@uva.nl 


Magnetars are highly magnetized neutron stars; their formation mechanism 
is unknown. Hot helium-rich stars with spectra dominated by emission lines 
are known as Wolf-Rayet stars. We observe the binary system HD 45166 using 
spectropolarimetry, finding that it contains a Wolf-Rayet star with a mass of 2 
solar masses and a magnetic field of 43 kilogauss. Stellar evolution calculations 
indicate that this component will explode as a type Ib or IIb supernova, and 


the strong magnetic field favors a magnetar remnant. We propose that the 


magnatized Wolf-Rayet star formed by the merger of two lower mass helium 


stars. 


Neutron stars form in supernovae, by the collapse of stellar cores that exceed the Chan- 
drasekhar mass limit [mass M = 1.4 solar masses ( M.;)]. Bare stellar cores can be exposed as 
hot, evolved helium-rich stars that have shed their outer hydrogen-rich layers. A subset of these 
massive helium stars are observed as Wolf-Rayet stars, which have spectra dominated by broad 
emission lines produced by strong stellar winds (/, 2). Massive helium stars (i.e. M = 1.4 Mọ) 
are thought to be stripped products of massive stars that lost their hydrogen-rich envelopes 
through stellar winds, eruptions, or interactions with a binary companion (3, 4). Alternatively, 
massive helium stars may be produced through the merger of lower-mass objects (5). 

Roughly 10% of young neutron stars have magnetic fields > 1014 gauss (G) (6). These 
are known as magnetars; their origin is debated (7, 8). One formation scenario invokes fossil 
magnetic fields rooted in the pre-collapse massive core (9). About 7-10% of massive main- 
sequence stars have strong (several kG) large-scale surface magnetic fields (/0, 1 1); these could 
be progenitors of magnetars. However, corresponding magnetic fields have not been detected 
in evolved massive stars (/2). Strongly magnetic low-mass helium stars have been observed 
(13-15), but not massive magnetic helium stars (i.e. exceeding the Chandrasekhar mass limit). 

The HD 45166 (ALS 8946) system is a binary comprising a main sequence star (classified as 
spectral type B7 V) with a hot stellar companion. The hot companion has a spectrum dominated 
by the characteristic emission lines of a Wolf-Rayet star (Fig. 1), but it was classified as "quasi 
Wolf-Rayet" (qWR) star owing to its peculiarly narrow emission lines (widths of hundreds of 
km s~' instead of thousands), spectral variability, and the anomalous presence of strong carbon, 
oxygen, and nitrogen lines in its spectrum. Previous radial velocity (RV) measurements have 
shown a 1.6d periodicity in the velocity of the B7 V component, interpreted as the orbital 


period of the system (76). This implied a mass of 4.2 + 0.7 Mọ for the Wolf-Rayet component 
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and a pole-on orbital configuration (inclination 7 = 0.7?) (76). That mass is well below the 
typical masses of Wolf-Rayet stars in our Galaxy [M = 8 Mo, (17)]; no other massive helium 
stars are known with M < 8 Mọ (18). Comparison to stellar atmosphere models that do not 
assume local thermodynamic equilibrium (non-LTE) found that the Wolf-Rayet component is a 
hot (surface effective temperature of T, = 70 kK), helium-rich star with enhanced nitrogen and 
carbon contents compared to a solar composition (79). A latitude-dependent wind model was 
invoked to reproduce the spectrum of the Wolf-Rayet component (179). The implied mass-loss 
rate is orders of magnitude higher than predictions for helium stars of this mass (20—22). 

The presence of a C and N emission-line complex in the range 4430 — 4460 À (referred to 
as the Of?p phenomenon) and strong spectral variability are typical signatures of magnetism 
in hot stars (25, 24). Here, we use spectropolarimetry to investigate whether the Wolf-Rayet 


component in HD 45166 is magnetic. 


Observations of HD 45166 


We collected eight high-resolution spectropolarimetric observations (spectra of Stokes param- 
eters J and V) of HD 45166 in February 2022 with the ESPaDOnS spectropolarimeter at the 
Canada-France-Hawaii Telescope (CFHT) (25). The spectra cover the range 3668 - 10480 À 
and are used to measure the magnetic-field strength in the Wolf-Rayet and B7 V components 
of HD 45166. We additionally use spectra acquired with three spectrographs for RV monitor- 
ing (25). 103 spectra were obtained with the Coudé spectrograph at the 1.6 m telescope of Lab- 
oratório Nacional de Astrofísica (LNA), mainly covering the range 4520 — 4960 À. 36 spectra 
were obtained with the Fiber-fed Extended Range Optical Spectrograph (FEROS) at the 1.52 m 
telescope of the European Southern Observatory (ESO), covering the range 3830 — 9215 Å. 
Finally, 28 spectra were acquired with the HERMES spectrograph mounted on the 1.2 m Mer- 


cator telescope, covering the range 3770-9000 A. The analysis of the spectrum and spectral 


energy distribution (SED) relies on additional ultraviolet (UV) and photometric data (25). 


Evidence for magnetism 


We detect strong circular polarization (Stokes V) in the ESPaDOnS spectra in the majority of 
lines associated with the Wolf-Rayet component (Fig. 2). We also detect the magnetic splitting 
(i.e. Zeeman splitting) of two spectral lines belonging to O v in our spectrum (Fig. 2C-D) that 
form in or close to the stellar surface. We measure the magnetic field from the separation of 
the split Zeeman components of these lines, finding a mean field modulus (B) wr = 43.0 + 
2.5 kG and a mean longitudinal field (B;)qwn = 13.5 + 2.5 kG (Table S1). The ratio of the 
longitudinal field to the modulus (1/3rd) is consistent with a dipolar magnetic field viewed 
near the magnetic pole (26). 

Such a strong magnetic field in the Wolf-Rayet component implies that its emission-line 
spectrum is formed in plasma confined to the magnetic field loops (the magnetosphere), and not 
in a radially expanding stellar wind (25). In that case, one-dimensional (1D) stellar atmosphere 
models [as used in previous studies (/9)] are insufficient for analyzing the emission lines, which 
affected the derived physical parameters. We reanalyze the spectra and SED using the non- 
LTE Potsdam Wolf-Rayet (POWR) model atmosphere code (27, 28). While the POWR code 
is 1D, we only use it to analyze features originating in the stellar surface of the Wolf-Rayet 
component, ignoring its emission features (25). The resulting physical parameters are listed in 
Table S1. We find a high effective temperature for the Wolf-Rayet component, though ~15kK 
lower than previously determined (79). The effective temperature (7T. = 56.0 + 5.0 KK) and 
bolometric luminosity (in units of solar luminosity Lo, log L/ Lc = 3.830 + 0.050) indicate 
that the magnetic Wolf-Rayet component is not a main sequence star, but an evolved object. We 
applied the same analysis to the B7 V component. By using its derived stellar parameters as an 


input in the BONNSAI Bayesian tool (29), we derive Mp = 3.38 +0.10 Mo and 105+ 35 Myr 
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for mass and age of the B7 V component, respectively (25). 

The Wolf-Rayet component exhibits changes in the strength of lines such as He II A4686. 
These appear to be periodic (Fig. 3A), with a best-fitting period of 124.8 d (Fig. 3B). Periodic 
changes in the line strengths of hot magnetic stars are typically interpreted as their rotational 
periods (30, 31), implying a rotational period Pj, qwn = 124.8 + 0.2d for the Wolf-Rayet 
component. This period is consistent with the projected rotational velocity indicated by the 


narrow O v lines (projected rotational velocity v sini < 10 km s71). 


The orbit of HD 45166 


The ESPaDOnS spectra indicate that the B7 V component exhibits spectral line profile varia- 
tions caused by non-radial gravity-mode pulsations (Fig. S9). This is supported by a frequency 
analysis of a light curve obtained with the Transiting Exoplanet Survey Satellite (TESS) per- 
formed in this study (25), which yields 1.6 d as one of the significant periods. we conclude 
that the 1.6-d period previously attributed to orbital motion (76) is caused by pulsations in the 
B7 V component. We therefore need to reassess the orbit of and the mass of the Wolf-Rayet 
component, M wr. 

To determine the orbit, we compiled spectroscopic data spanning 24 years taken with the 
facilities listed above (25). These data show a long-term anti-phase motion of the B7 V and 
Wolf-Rayet components (Fig. 4). We find multiple periodicities associated with both the Wolf- 
Rayet and B7 V components (25), but we determine a best-fitting orbit (Table S1). We find an 
orbital period P = 8200 + 190d and a semi-major axis a = 10.5 + 1.8 astronomical units 
(au). This indicates that the components are much more widely separated than had previously 
been determined. Our derived mass ratio of q = Mqwn/Mg = 0.60 + 0.13, combined with 
the mass of the B7 V component (Mg = 3.38 + 0.10 Mo) we derived above, implies Mawr = 


2.03 + 0.44 Mo (uncertainties are 68% confidence intervals). This is less than the 4.2 Mo 


previously reported, but still above the Chandrasekhar limit. The derived mass and luminosity 
are consistent with mass-luminosity relations for He stars (32). The large separation and implied 
orbital inclination of i = 49 + 11? (derived from My sin? i and Mp, Table S1) explain the low- 


amplitude RV motion of both binary components, without invoking a pole-on configuration. 


Implications on the magnetar population 


With a mass of 2.03 + 0.44 Mo, we can expect the Wolf-Rayet component to evolve until it 
collapses into a neutron star. This is supported by evolution models constructed for the system, 
although the final fate of the Wolf-Rayet component depends on modeling uncertainties (25). 
Upon core-collapse, magnetic flux conservation leads to an increase in the magnetic field at the 
surface. With a stellar radius of R, qwg = 0.88 +0.16 solar radii (Rọ) calculated via the Stefan- 
Boltzmann relation, (B)qwn = 43.02: 2.5 KG, and assuming a final neutron-star radius of 12 km 
(33), we calculate a final magnetic field of the neutron star (B) jg = (1.11 + 0.42) x 101G. 
This is within the range observed for magnetars [(B) => 10!^ G, (34)]. Our measurements 
and evolution models therefore indicate that the Wolf-Rayet component in is an immediate 
progenitor of a magnetar. 

All magnetars in the Milky Way are isolated i.e. do not have a binary companion (6). For 
the Wolf-Rayet component in HD 45166, we expect the mass-loss and velocity kick imparted 
on the magnetar by the supernova explosion to disrupt the system, given the large orbital sepa- 
ration. With an estimated rotation period of 125 d and an estimated radius of ~ 0.3 Re for the 
helium core of the Wolf-Rayet component (35), angular momentum conservation implies that 
the magnetar immediately after collapse would have a spin period < 40 ms. This is similar to 
the spin period of the famous Crab Pulsar (33 ms), which is a neutron star formed about 1000 
years ago. Such a spin is not expected to provide sufficient energy to power a superluminous 


supernova or a long duration gamma-ray burst (36, 37). 


Evolutionary history of the system 


We next consider how the Wolf-Rayet component itself formed. We exclude the possibly that 
it is the stripped descendant of a massive star, because to produce a 2 Mo stripped core, the 
progenitor star would need to have had an initial mass of ~ 10 Mo. Single-star evolution models 
do not predict that stars of that mass strip themselves (22, 35), and the B7 V companion is too 
far away for binary interactions to have stripped it. In addition, the total lifetime (including post 
main-sequence evolution) of a 10 M; star would be ~ 30 Myr (39), well below the derived age 
of the B7 V component (105 + 35 Myr), so we reject the possibility that they could both be 
present in the same binary system. 

Stellar mergers have been proposed as a potential origin of magnetic stars (40—43). Strong 
magnetic fields have been identified in low-mass helium stars (OB-type subdwarfs) and were 
suggested to originate in merger events between two white dwarfs (13, 14). However, to produce 
a ~ 2 Mo helium star through a merger of white dwarfs, the merger would need to have involved 
rare, massive CO or ONe white dwarfs. Models predict that such a merger product would to 
either immediately explode as a supernova (44), or do so after a short lifespan of approximately 
10 kyr (45), suggesting that observing a stellar product following the merger is unlikely. 

We therefore propose that the Wolf-Rayet component in formed from the coalescence of the 
helium cores of two intermediate-mass stars that were bound in a close binary. We construct an 
evolutionary model of this scenario using the Modules for Experiments in Stellar Astrophysics 
(MESA) stellar evolution code (46) (Fig. 5). The system is started in a triple configuration with 
a tight inner binary and a distant tertiary representing the B7 V component. We find that the 
primary (more massive) star of the inner binary expanded and interacted with its companion, 
thereby losing its outer layers and becoming a stripped star, while the companion accreted 


material and became rejuvenated with hydrogen. Later, the secondary star expanded, leading to 


an unstable mass-transfer. This process leads to the formation of a gaseous envelope around the 
two stars, during which the two stars lose orbital energy due to friction with the envelope and 
spiral inwards (common-envelope evolution). Due to the high binding energy of the hydrogen 
rich layers, this phase ends with both helium cores merging into a ~ 2 Mo magnetized helium 
star while most (but not all) of the hydrogen envelope is ejected. Being roughly 200 times more 
luminous than a main sequence 2 Mo star (47), the merger product has a high luminosity-to- 
mass ratio. This, combined with its high effective temperature, launches a radiatively-driven 
outflow. In the absence of a magnetic field, such an outflow would not be easily detectable in 
a spectrum. However, the trapping of the outflow by the magnetic field is known to increase 
the density of the circumstellar material and produce spectral emission lines like those observed 
in the Wolf-Rayet component (26, 48). The distant tertiary star in our model does not impact 
the final outcome of the inner binary, except for perhaps catalyzing the merging of its stellar 
components (49). Given its large orbital separation, we do not expect the tertiary star to show 
evidence for accreted material from the merger ejecta. 

The scenario proposed here is quantitatively and qualitatively consistent with the observed 
properties of the system. The MESA model reproduces the masses of the two components and 
the age of the system. The merger provides an explanation for the emergence of a magnetic field 
in the Wolf-Rayet component, and the magnetic field provides an explanation for the presence 
of emission lines in the spectrum of a 2 Mọ helium star. The binary nature of enabled us to 
use the companion as a clock to constrain the evolutionary path of the system, and as a scale to 
determine the mass of the Wolf-Rayet component — conditions that are rarely present in other 
proposed merger products (42, 50). Given the proximity of HD 45166 to Earth (~ 1 kiloparsec), 
other massive magnetic helium stars have likely already been spectroscopically identified as 


Wolf-Rayet stars but not recognized to be magnetic (see Supplementary Text). 
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— — 8 Dec 2020 (MJD 59191.2) 
---- 13 Feb 2021 (MJD 59258.1) 


Of?p phenomenon 
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Figure 1: Indications for magnetism in HD 45166. Two optical HERMES spectra of 
HD 45166 were taken ~ 70d apart (see legend). The presence of the so-called Of?p phe- 
nomenon (23, 24) and the strong variability observed in the He I1 \4686 line are common fea- 
tures of hot magnetic stars, indicating that the Wolf-Rayet component might be magnetic. The 
narrow Mg II A4481 line is due to the B7 V component. 
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Figure 2: Stokes spectra of HD 45166. The panels show the intensity spectrum (7, black lines), 
diagnostic null spectrum (N; blue lines), and Stokes V spectrum (upper curves, red lines) of the 
co-added ESPaDOnS spectrum. The N and V spectra have been vertically shifted and scaled 
for display purposes (the multiplication factors are listed in each panel). (A) — (H) Several 
diagnostic lines of the Wolf-Rayet component. (I) An O I triplet associated with the B7 V star. 
Zeeman splitting is visible in the O v A4930 and O v \5114 lines (C, D). We use these lines to 
measure the magnetic field strength. There is no Stokes V signature visible for lines associated 
with the B7 V star (panel I). 
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Figure 3: Variability of the He I1 \4686 emission line. 
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(A) HERMES spectra of the 


He II A4686 line (color depicts the time of observation), illustrating the changing strength of 
this line on a timescale of weeks. (B) Line strengths (equivalent widths) of He II A4686 over 
the 24-yr spectroscopic dataset. Different colors and plotting symbols indicate different spec- 
trographs, as indicated in the legend. The measurements were phased at the derived period of 
124.8 d. We interpret this as the rotational period of the Wolf-Rayet component. Error bars are 


estimated lo confidence intervals. 
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Figure 4: Two-component orbital solution of HD 45166. The individual RVs of the Wolf- 
Rayet and B7 V components (gray symbols) are binned by 30d (green and black symbols 
for Wolf-Rayet and B7 V components, respectively). The different instruments are shown as 
distinct symbols (see legend). Solid lines correspond to the best-fitting RV curves for the Wolf- 
Rayet (black) and B7 V (green) components. The fit yields a reduced x? of 1.21. Error bars 
indicate estimated lo confidence intervals (25). 
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Figure 5: Models for our proposed evolutionary scenario for HD 45166. MESA models are 
plotted on a T, — log L diagram for a triple system comprising an inner binary with an initial 
period of 2 d and initial masses of 5 M; (primary, blue line) and 3 Mọ (secondary, orange line), 
and a tertiary (gray line) with an initial mass of 3.4 Mo. Solid lines correspond to the evolution 
of the three components before the inner binary merges; dashed lines represent the post-merger 
evolution. Dots along the evolution tracks are separated by 0.1 Myr of evolution. The merger 
position is marked with a red cross. Yellow and red highlighted regions correspond to mass 
transfer on the main sequence (case A) and after the main sequence (case AB), respectively. 
The upper brown and pink tracks depict the evolution of the merger product (representing the 
Wolf-Rayet component) for common-envelope ejection efficiencies of acg = 0.25 and 0.5, 
respectively (25). Black plus signs mark the observed locations of the Wolf-Rayet and B7 V 
components. Purple filled circles mark the core He depletion of the merger, which occurs after 
133.2 and 133.7 Myr for acg = 0.25 and 0.5, respectively. 
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Materials and methods 
Observational data 


We acquired spectropolarimetric observations (Stokes / and V) of HD 45166 on 2022 February 
17 and 19 with the ESPaDOnS spectropolarimeter (53) mounted on the CFHT (program ID 
22AC08; PI: Shenar). A total of eight 1-hr exposures (four per night), each divided into four 
sub-exposures, resulted in typical signal-to-noise (S/N) of ~ 150 — 200 per wavelength bin at 
~ 5000 À and a resolving power R = 65000. The spectra cover the wavelength range 3668 
— 10480 Å. The data were reduced and wavelength-calibrated using CFHT's standard UPENA 
pipeline (54). Including the subexposures, a total of 32 Stokes J and V ESPaDOnS spectra 
were obtained. The co-added ESPaDOnS / and V spectra are shown in Fig. 2. 

Additionally, we measured RVs from archival optical spectra obtained with three additional 
spectrographs. The first set of archival data was acquired between 1998 to 2004 using the Coudé 
spectrograph at the 1.6 m telescope of the LNA in Itajubá, Brazil (55), and has been previously 
described by (76). The majority of LNA spectra cover a modest spectral range (4520 — 4960 A) 
at a resolving power of R ~ 7000, reaching a typical S/N of 20 — 30 in the continuum. A few 
additional higher-resolution spectra (R ~ 20 000) covering the range 3830 — 4580 Å were also 
used. In total, 103 LNA spectra were available. 

The second set of archival spectroscopic data were obtained in 2002 with the FEROS in- 
strument (56) mounted at the 1.52 m telescope of the European Southern Observatory (ESO) in 
La Silla, Chile (56). These data cover a spectral range 3830 — 9215 À at a resolving power of 
R = 48000, with typical S/N of 80 in the range 4500 — 7000 A decreasing down to z 35 at the 
blue and red edges. A total of 36 FEROS spectra are available, taken at a high cadence over a 
few nights. The data reduction and calibration were also described previously (/6). 


The third and final archival spectroscopic dataset was acquired in 2019 — 2021, with the 
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exception of a single spectrum in 2012, using the HERMES spectrograph (57) mounted on 
the 1.2 m Mercator telescope at the Observatorio del Roque de Los Muchachos on La Palma, 
Spain (57). HERMES spectra cover the wavelength range from 3770 to 9000 À with a spectral 
resolving power of R 2:85 000. A total of 28 HERMES spectra were available, having a typ- 
ical S/N of 30 in the continuum. Standard reductions including bias and flat-field corrections 
and wavelength calibrations were performed using the HERMES data reduction pipeline (57). 
Barycentric corrections were applied. 

All spectra were normalized by fitting a piecewise linear function through identified con- 
tinuum regions selected homogeneously for all datasets. Given the low-amplitude motion of 
the binary components, stable wavelength calibration is essential. Relative wavelength calibra- 
tion was brought to a precision of ~ 0.1 kms! by cross-correlating the interstellar sodium 
lines Na 1 AA5890, 5896 in the FEROS, HERMES, and ESPaDOnS spectra with a selected ES- 
PaDOnS spectrum (which is wavelength-calibrated). However, this was not possible for the 
LNA spectra, which do not cover these lines. We therefore consider the RVs of the LNA data 
as less reliable, and allow for a systematic RV shift between them and the other datasets (see 
below). 

For the spectral analysis, we append the optical spectra with a set of UV spectra retrieved 
from the Mikulski Archive for Space Telescopes (MAST). The data were acquired with the 
International Ultraviolet Explorer (IUE) (program IDs: WRJSH, JA017; PIs: Heap, Willis, 
respectively) and presented previously (58, 59). The spectra cover 1150 — 2150 À at a resolving 
power of R œ~ 10000. There are substantial changes in P-Cygni lines that occur over short 
and long timescales (58, 59). We use a co-added spectrum formed from 20 high-cadence IUE 
spectra obtained consecutively on 1988 February 3, 4, and 5 for 36hr. The individual spectra 
have a S/N of 5-6, such that the co-added spectrum has a total S/N of « 20. While the process 


of co-adding smears the small-amplitude variability of the Wolf-Rayet component on the 5 hr 
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period (see below), the impact is negligible compared to the spectral resolution. Since we 
only use this spectrum to estimate the effective temperature of the Wolf-Rayet component and 
analyze the SED, the process of co-adding the data has no impact on our results. The spectrum 
is normalized using the model continuum obtained from the SED analysis (see below). 

For the analysis of the SED, we also retrieved from MAST a far-UV spectrum covering the 
spectral range 900 — 1200 À obtained with the Far Ultraviolet Spectroscopic Explorer (FUSE), 
which have also been previously described (60). The spectrum was acquired on 2001 Mar 5 (PI: 
Willis, program ID: P224), has a resolving power of R z 20000, and a S/N of ~ 100. Since the 
FUSE and IUE spectra overlap, we use the IUE spectra to calibrate the FUSE spectrum rather 
than the standard flux calibration, such that the overlapping regions (~ 1150 — 1200 A) of the 
IUE spectra and FUSE spectra agree, resulting in a multiplicative factor of 0.8 of the FUSE 
spectrum. For the SED analysis, we append the FUSE and IUE spectra with a low-resolution 
flux calibrated IUE spectrum obtained on 1982 Sep 06 (PI: Stickland, ID: EI273), which covers 
the range 1850 — 3350 A at a resolving power of R ~ 500 and S/N% 5. 

Finally, for the SED analysis, we use U photometry (61), BV photometry (62), JH K pho- 
tometry (63), and Wide-field Infrared Survey Explorer (WISE) photometry (64). 


Analysis of the spectropolarimetric data 


Stokes V variations are detected in all ESPaDOnS spectra in multiple lines attributed to the 
Wolf-Rayet component (Fig. 2), most prominently lines of He II and lines of highly ionized 
elements. Most lines are in emission in the intensity spectrum, although O v 45114 appears in 
absorption. No similar Stokes V features are observed in the narrow absorption lines belong- 
ing to the B7 V component (Fig. 2I). Since the four individual ESPaDOnS Stokes V spectra 
obtained during the two observing nights do not show notable variability in the context of the 


magnetic properties, they were co-added to achieve a higher S/N. 
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The circular polarization features are qualitatively consistent with the Zeeman effect due to 
the presence of a strong magnetic field in HD 45166. While the features exhibit some morpho- 
logical diversity (Fig. 2), they are all - except for O v A4930 - of the same polarity, with positive 
circular polarization observed in the blue wing, and negative polarization observed in the red 
wing. O v A4930 is a weak emission line that appears distinct from the other emission features, 
with a polarity opposite to that of all other lines in the spectrum (Fig. 2D). 

The O v 45114 and 4930 lines indicate the magnetic nature of the Wolf-Rayet component 
and provide quantitative information about its magnetic field. Both lines exhibit resolved Zee- 
man splitting into triplets that correspond in wavelength to the locations of the extrema and 
inflection points of their associated Stokes V features. Their narrow profiles, and the pure ab- 
sorption profile of the O v 45114 line, indicate that they form in or close to the stellar surface. 
The inverted Stokes V signature of O v \4930 can be explained by the Zeeman effect operating 
on a spectral line formed in emission, and not in absorption. As discussed previously (65) for 
the strongly magnetic O-type star NGC 1624-2, spectral lines of hot stars can form in emission 
under non-LTE conditions even in or close to the stellar surface. 

To quantify the strength and geometry of the detected magnetic field, we measured the mean 
magnetic field modulus (5) and mean longitudinal magnetic field (B,) from the O v AA5114 
and 4930 lines. We renormalized the two lines to the local continuum. We measured the mean 
field modulus by fitting Gaussians to the three components of each Zeeman triplet to measure 
the average separation between the two polarized line components and the unpolarized line, 
AXzeeman. We then compute (B) = AXzeeman/(4.67 x 107 g AS), where g) = 1.0 is the 
theoretical Landé factor of the line and \y = 4930 À or 5114 À is the rest wavelength. We 
find (B)qwn = 44.1 + 2.0kG and (B)qwn = 41.8 + 1.0 kG for the O v \ 4930 and A5114 
lines, respectively. The shallow depths of the two O v lines and uncertainty in the position- 


ing of the continuum leads to systematic errors that are not included in the formal lo uncer- 
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atinties provided, which we thus ignored when computing the mean and standard deviation: 
(B)qwn = 43 + 2.5 KG. We measured a mean longitudinal field of (B;)qwn = 16.1 + 2.0kG 
and (B;)qwn = 11.1 + 1.0 kG for the O v A4930 and A5114 lines, respectively, using the first- 
order moment method (66), adopting the same Landé factor and rest wavelengths. Taking the 
mean and standard deviation yields (B;)qwg = 13.5 + 2.5 KG. As before, we do not weigh the 
values by the errors due to unaccounted systematic errors. 

We also extracted Least-Squares Deconvolution (LSD) mean profiles of the B7V companion 
(Fig. S1). We employed a standard solar abundance B7 line mask (/0) and used the ILSD 
code (67, 68). From the LSD profile we determined the detection probability in both Stokes V 
and the diagnostic null (both non-detections), and determined the upper limit of the star's mean 
longitudinal magnetic field [< 60 G at 3c confidence; (66)]. 

We perform a more general comparison between the Stokes / and V spectra of the Wolf- 
Rayet component and NGC 1624-2 (Fig. 52). NGC 1624-2 is a magnetic star showing the 
Of?p phenomenon (23, 24) with a very strong magnetic field of ~ 20kG (65), sufficiently 
large to produce marginal Zeeman splitting of some Stokes / and Stokes V profiles of some 
individual spectral lines. The star is known to have a large, dense magnetosphere (69, 70). The 
lower effective temperature of NGC 1624-2 of 35 + 2 kK (65) leads to absorption lines of He I 
in its spectrum. As a consequence, essentially all H and He I lines, whether in emission or 
absorption, show Zeeman signatures. In HD 45166, however, no Zeeman signature associated 
with He I lines is seen, even though He I lines are present (in emission) in its spectrum (Fig. 2H). 
Instead, Zeeman signatures are detected in association with numerous lines of highly ionized 
light elements, in particular C IV, N IV, N V, and O v. Another difference is in the ratio of 
the He I1 A4686 line to the Ha + He 11 A6560 blend. While the Ha + He 11 46560 blend is of 
similar strength in both stars, the He I1 A4686 line is almost 20 times stronger in the Wolf-Rayet 


component of HD 45166, indicating that it is helium rich. 
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The observed properties of the magnetic field in the Wolf-Rayet component are similar to 
those found in other hot magnetic stars with dipolar magnetic fields (26), indicating that the 
Wolf-Rayet component of also hosts a global magnetic field. In this configuration, the polar- 
ized flux is largely generated near the stellar surface, where the magnetic field is strongest. As 
in other magnetic star, the emission lines then originate from stellar wind trapped by the closed 
magnetic loops (the magnetosphere). The Stokes-V Zeeman signature is produced predomi- 
nently in the surface layers of the star. This is because the Zeeman signature is seen only in 
high ionization lines, which are formed in or close to the stellar surface according to the model 
atmospheres computed for the star (see below). At larger radial extent, no or weak Zeeman 
signature is seen, both because the magnetic field strength rapidly drops radially (x r^? in 
a dipolar field, r being the radial distance), and because the signature is diluted in the larger 
emitting volume, which exhibits opposing polarities. While some contribution to the observed 
circular polarization could be produced by field lines extending into the star's wind, the lack of 
any circular polarization detected in the emission lines of He I indicates that such contributions 
are very small. In Fig. S3 we show LTE polarized spectrum synthesis calculations with the 
ZEEMAN code (51, 71, 72) for the O v \5114 line, adopting a dipolar magnetic field of 45 kG 
with the magnetic field axis parallel to Earth's line-of-sight. The narrow Zeeman components 
constrain the non-thermal line broadening to 10 km/s. In Fig. S3 we have adopted a Doppler 
width of 6 km/s to account for rotational and turbulence broadening. Although we have adopted 
a simplified magnetic geometry, the Stokes J and V profiles are reproduced by this basic model. 

The Stokes V spectrum thus indicates that the underlying intensity spectrum of the Wolf- 
Rayet component, which is largely hidden behind the magnetosphere, is characterized by the 
presence of lines belonging to highly ionized species such as N V, O v and the lack of lower- 
ionization lines such as He I and O 11. To fulfill such conditions, the Wolf-Rayet component 


needs to be hot (T, = 45 kK), consistent with our spectral analysis (see below). 


a 
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Spectral analysis 


To derive the stellar parameters of the components, we use the non-LTE Potsdam Wolf-Rayet 
(POWR) model atmosphere code (27, 28, 52, 73), which performs radiative transfer in an expand- 
ing atmosphere. Models are defined by the chemical abundances, the surface effective temper- 
ature 7,, the bolometric luminosity L, the surface gravity g, and wind parameters such as the 
mass-loss rate M and terminal wind speed væ. The radius R, is given by the Stefan-Boltzmann 
relation, R, x VLT ig In POWR, the effective temperature is defined with respect to a contin- 
uum Rosseland optical depth of TRoss = 20 as a proxy for the hydrostatic layers. These layers 
are deeper than the photosphere (Tross = 2/3), at which the photospheric effective temperature 
Tog is defined. For a star with a negligible wind (such as the B7 V component), T; ~ Teg. For 
the Wolf-Rayet component, however, Ter is the effective temperature with respect to a layer 
that is embedded in the circumstellar material. Hence, we use 7’, when comparing to evolution 
models. 

In the subsonic regime, the density is determined by hydrostatic equilibrium. In the super- 
sonic regime, the radial velocity field is assumed to follow a 5-law (74), and the density follows 
from the continuity equation. The supersonic domain does not impact our analysis, since the 
wind of the B7 V component is negligible, and since the emission features of the Wolf-Rayet 
component are ignored. In the main non-LTE iteration, the line profiles are assumed to have 
a constant Doppler width of 30 kms" !. During the integration, the Doppler widths are deter- 


mined from the depth-dependent thermal broadening and microturbulence € (75). 
Analysis of the B7 V component 


Atmospheric abundances of the B7 V component are fixed to solar (76), which is consistent with 
the data. The surface gravity of the B7 V component cannot be determined, since the diagnostic 


Balmer lines are entangled with magnetospheric emission and absorption originating in the 
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Wolf-Rayet component in a non-trivial manner. We therefore assume log g/ cms ? = 4.0 (47); 
varying log g by 0.2 dex has little impact on the remaining diagnostic lines. Wind parameters 
are set to negligible values given the lack of wind diagnostics. Microturbulence is set to £ — 
5 kms‘. We determine T, and the overall light contribution of the B7 V component in the V- 
band, [p (V), using diagnostic lines of O 1, Mg II, and Si 11, avoiding He I due to possible cross- 
contamination by the magnetosphere of the Wolf-Rayet component (Fig. 54). The projected 
rotational velocity v sin ? and macroturbulent velocity Umac are derived using a Fourier method 
(77), implemented by the IACOB-BROAD tool (78—80). However, in the domain where v sin? 
is small, this method is known to lead to degeneracies (79), and we therefore only provide 
upper limits in Table S1. The rotation and macroturbulence are accounted for by convolving the 
synthetic spectrum of the B7 V component with corresponding rotational profiles with the value 
corresponding to the upper limit given in Table S1. the projected rotational velocity of the B7 V 
component (vsini < 10 kms?) is low for this spectral type (87), though other comparably 
slow rotators exist. While a pole-on configuration of the rotational axis of the B7 V component 
would provide an explanation for this, such a configuration would not be consistent with the 
orbital inclination derived below (i = 49 + 11°), unless it the orbital and rotational axes are 


misaligned. More likely, the B7 V star is intrinsically a slow rotator. 
Analysis of the Wolf-Rayet component 


For the Wolf-Rayet component, we adopt previously derived elemental abundances for He, C, 
N, and O (19); Xy, = 0.67, Xc = 0.0059, Xn = 0.002, and Xo = 0.0015, which leaves 
Xy = 0.32 (values are mass fractions of the respective element). The remaining abundances 
are set to solar. Implications of this abundance pattern have been discussed previously (/9); it is 
anomalous compared to either Wolf-Rayet stars or central stars of planetary nebulae. However, 


these abundances were derived assuming that the emission features stem from a radially expand- 
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ing outflow, disregarding the magnetic confinement of the wind. We therefore refrain from inter- 
preting the abundances quantitatively. Unlike the previous investigation (79), we rely solely on 
the iron forest in the UV, which form in the stellar surface and are not contaminated by the mag- 
netosphere. For this purpose, the co-added IUE spectrum is used, which is normalized using the 
added model continua of both stars. We adopt a mass-loss rate of log M /( Mz yr!) = —9.0 for 
the Wolf-Rayet component; The high mass-loss rate of log MK Mz yr~') = —6.7 previously 
derived (79) was a consequence of ignoring the magnetosphere, so we do not adopt it. We use 
regions dominated by lines belonging to Fe Iv, v, VI (Fig. S5). By using the UV lines, the con- 
tribution of the B-type companion is negligible (see SED analysis below). The narrow absorp- 
tion features in the visual imply a combined projected rotational velocity and microturbulence 
vsini,£ < 10 kms" !, and we adopt these upper limits for the computation of the spectrum. 
The low projected rotational velocity of the Wolf-Rayet component can be understood as a con- 
sequence of angular momentum loss via the magnetic field (82, 53). From the iron forest, we 
determine T, = 56 + 5 KK, which is lower than the value of 70 kK derived previously (79) by 


analyzing magnetospheric features using a wind model. 
Spectral energy distribution and age 


Using the POWR models of the two components and their light ratio in the visual, we compute 
the luminosities and reddening by fitting the multiwavelength SED to the sum of the models 
(Fig. S6). We use a reddening law (84) with a total-to-selective extinction ratio of Ry = 3.1. 
We find log Lqwn/ Lc = 3.830 + 0.050 and log Lp/ Lco = 2.250 + 0.050 [Lo], with a color 
excess of Ep. y = 0.210 + 0.010, which is consistent with the strength of the interstellar lines 
(e.g., Lyo). The parameters derived for the B-type component are consistent with its spectral 
type (47). A flux excess is noted in the infrared which grows gradually between 2 um and 


22 um; it could originate in free-free emission stemming from the trapped outflow, which is not 
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included in our model. 

We used the BONNSAI Bayesian tool (29, 85) to estimate the current mass and age of 
the B-type component using single-star evolution tracks (39), using our derived T,, log L, and 
v sin i as input parameters. We find a mass of Mp = 3.388+0.10 Mo and an age of 1053-35 Myr, 


which we consider as the age of the binary. 


TESS light curve and the 1.6 d period 


To analyze the light curve of HD 45166, we downloaded 15 x 15 pixel cutouts from the TESS 
full frame images (FFIs) stored at MAST using the ASTROCUT package (56). We extracted 
light curves for HD 45166 by defining custom aperture masks, and estimated the background 
flux from the median flux per frame excluding pixels that contain stellar flux. We subtracted 
the background flux from the target flux, normalized by dividing through the median. A final 
linear trend was subtracted from each light curve to remove remaining instrumental effects. Our 
FFL-extracted light curves for sectors 6 and 33 are shown in the left panel of Fig. S7A and B. 
To determine the variability in each light curve, we used Fourier analysis and the method of 
iterative pre-whitening (87). We fit the light curve of each sector with a model function of the 


form: 


Am (t) — ` Aj cos (27 vj (t — tree) + ài), (S1) 


where t is the time stamps of the light curve with units of BJD—2457000.0, with BJD being 
the Barycentric Julian date (preferred over MJD for high-precision space-mission photometry), 
tref is a reference time at the end of the light curve of each sector (trer = 1490.0 and 2225.0 for 
sectors 16 and 33 in units of BJD-2457000.0, respectively), and v;, A;, and ¢; are the frequency, 
amplitude, and phase of the of the itè cosinusoid term, respectively. 


To determine the total number of significant frequencies, we use a significance criterion of 


33 


an amplitude signal-to-noise ratio (S/N) larger than four in the amplitude spectrum, where the 
signal is defined as the amplitude of a peak and the noise is defined as the average amplitude 
in the residual amplitude spectrum in a symmetric window of width 1 d^! centered at the ex- 
tracted peak for each iteration (87). In total, we detected two and three frequencies fulfilling 
our significance criterion in the sectors 6 and 33 TESS data of HD 45166, respectively. We 
performed an independent multi-frequency non-linear least-squares fit to each light curve us- 
ing Eq. (S1) to optimize the frequencies and their corresponding amplitudes and phases, and 
determine correlated uncertainties, which are provided in Table S2. 

The shorter 10-min cadence of the sector 33 data yields a lower noise level in the amplitude 
spectrum, on average, because of the factor three larger number of data points compared to the 
sector 6 data. This contributes to the detection of the additional frequency in the sector 33 TESS 
data, and the overall smaller uncertainties for the optimized parameters in Table S2. Given the 
low frequency resolution of a short (27 d) TESS light curve of 1/A(t) ~ 0.04 d^!, the two 
lower frequencies in sector 6 are the same as those in sector 33. 

Multi-periodic variability in a period range of several days is typical of gravity mode pulsa- 
tions in main-sequence B-type stars (88). The dominant frequency of ~ 0.6 d^! lies close to the 
period of 1.66 d detected from RV variations for the B-type component [ (/6), see also below]. 
We infer that the 1.6 d period, previously attributed to the orbital period of the system, is instead 
due to pulsations in the B7 V component. This is consistent with the much longer orbital period 


we find below. 


Revised orbital analysis of HD 45166 


To measure the RVs from the available spectra, we cross-correlate the observations with ap- 
propriate templates, following an established procedure (89). Measuring the RVs of the B-type 


component is only possible with the FEROS, HERMES, and ESPaDOnS spectra given the lower 
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wavelength coverage and quality of the LNA spectra. To form the template of the B7 V compo- 
nent, we first use the synthetic POWR model (Fig. 54). We then shift-and-add all observations 
using the measured RVs, and use this to produce a calibrated template of co-added observa- 
tions (90). The new template is then cross-correlated with the observations to yield the final 
absolute RVs of the B7 V component. We use the strong absorption lines Sill \\4128, 4131 
and Mg II AMASI.1, 4481.3 for the cross-correlation, and find that our conclusions do not de- 
pend on the choice of lines. 

For the Wolf-Rayet component, we use a similar procedure. However, lacking a synthetic 
template, we use one of the FEROS observations as a first template. We use the emission-line 
complex in the range 4630— 4660 À for cross-correlation. After measuring the preliminary RVs, 
which are measured relative to the observation used as a template, we form a co-added, high 
S/N template to remeasure the RVs. We allow for an overall systematic shift between the RVs 
of the Wolf-Rayet and B7 V components, such that they result in the same systemic velocity 
(see below). The choice of emission lines used for cross-correlation does not affect our results, 
though usage of highly variable lines such as He 11 A4686 yields a substantially larger scatter, 
and is therefore avoided. A total of 199 and 96 RVs for the Wolf-Rayet and B7 V components 
are available. An overview of the derived RVs is shown in Fig. S8, and they are also provied as 
supplementary data (table S1). 

Both components exhibit substantial variability on several timescales, and this variability 
appears to be associated with several periods, consistent with previous investigations (/6). The 
RVs of the B7 V component vary on a short timescale, and we can confirm the 1.6 d period 
found previously (76) and formerly attributed to the orbital period of the system. By subtract- 
ing several periods from the RVs measured for the Wolf-Rayet component, the 1.6 d period was 
previously found (76) in the Wolf-Rayet component as well. However, having repeated this ex- 


periment using our dataset, we do not confirm this result. Instead, we find that the 1.6 d period 
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is associated with non-radial gravity-mode pulsations in the B7 V component (see above). This 
is evident in the ESPaDOnS spectra, which have pronounced variability between the two ob- 
serving nights. To boost the S/N, we co-added all observations acquired during a single night, 
and show the resulting spectra in Fig. S9. This spectral variability is typical for non-radial 
gravity-mode pulsations in slowly rotating stars (97). 

The strong emission lines of the Wolf-Rayet star display variability with a short period of 
the order of 5 hr, as previously identified (76). The origin of this period, which persist pseudo- 
periodically for at least 20 yr, is not clear. However, this period cannot correspond to orbital 
motion, since the associated RV semi-amplitude is ~ 1 to 2 kms! and varies from epoch 
to epoch. The combination of a 5 hr period with the low RV amplitude would either imply a 
highly unlikely pole-on geometry, or a sub-stellar mass for the companion, which is incosinstent 
with stellar evolution. It is possible that the 5 hr period is due to pulsations in the Wolf-Rayet 
component, in a phenomenon related to its magnetic field, or in an interplay between the two. 

Upon inspection of the long-term variability of the RVs of both components, anti-phase mo- 
tion between the two components is apparent. While the short-term variability partly masks this 
trend, we conclude that the two stars are bound on a long-period orbit. In the HERMES and 
ESPaDOnS observations taken during 2019-2022, the RVs of the B7 V component in Fig. S8 
are persistently red-shifted with respect to those of the Wolf-Rayet component. However, ex- 
actly the opposite is seen in the FEROS spectra, acquired in 2002. The individual HERMES 
observation taken in 2012 again shows a reversal of the RVs. 

To obtain an orbital solution, we binned the RVs on a 30d baseline to smooth-out short- 
term variability in the two components. To account for the short-term variability, we attribute 
an additional uncertainty of 3 kms ! to the individual measurements (representative of the 
amplitude of the short-term RV variation), but divide this by the number of measurements per 


30d bin, though we require the uncertainties to remain at least 1 km s-1. We fitted the binned 
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RVs with a model of the RV curve resulting from the orbital parameters (P, to, e, w, K1, K2, vo) 
using the Python Imfit package (92, 93). We also included a systematic shift between the RVs 
of the Wolf-Rayet component and a systematic shift of the LNA RVs (see above). The LNA 
offset is found to be —2.8 + 1.6 kms !, which is small and does not impact our solution. The 
solution is shown in Fig. 4, which corresponds to the parameters given in Table S1. The derived 
period is 22.5 yr, though this value depends on the treatment of the uncertainties. The formal 
uncertainties on the orbital parameters are derived from the covariance matrix; a Markov Chain 
Monte Carlo (MCMC) simulation yielded similar uncertainties. 

While the ratio K2/K, yields q = 0.58 + 0.16, we can determine the mass ratio inde- 
pendently of the orbital elements by considering that, at each point of time, the RVs of the 


components scale in proportion to their semi-amplitudes (94): 


RV(B) = -q- RV(qWR) + C, (S2) 


where RV(B) and RV(qWR) are the RVs of the B7 V and Wolf-Rayet components at a given 
epoch, respectively, q = Mqwn/Mg is the mass ratio, and C is a constant that depends on the 
systemic velocities and q. By performing a linear regression to the measurements (Fig. S10), 
we find q = 0.60 + 0.13, which agrees well with the value derived from the orbital analysis. 
The evolutionary mass of the B7 V component, in combination with the mass ratio, leads 
to the mass of the Wolf-Rayet component. Using our derived evolutionary mass of Mg = 
3.38 + 0.10 and the derived mass ratio of q = 0.60 + 0.13 yields Mawr = 2.03 + 0.44 Mo. 
From Mg sin? i = 1.50 + 0.74 Mo and Mg = 3.38 + 0.10, we infer i = 49 + 11°. The long 
period of the system provides an explanation for the low RV amplitudes observed in the system, 
without the need to invoke a pole-on configuration. We find that the Wolf-Rayet component 


exceeds the Chandrasekhar mass limit, and is expected to undergo core-collapse (see below). 
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Evolutionary scenario 


Assuming that the Wolf-Rayet component formed through a merger of two helium stars (see 
main text), we construct an evolutionary model for the system using the MESA stellar evolution 
code (46). We considered the evolution of a 5M, star with a 3M, companion at an initial 
orbital period of two days, representing the inner binary in an initial triple. The more massive 
component overflows its Roche lobe first, undergoing mass transfer both during and after its 
main sequence, resulting in a stripped 0.63Mo star. Mass transfer increases the orbital period 
and separation to 66 days and 138 R5, while the accretor mass increases to 7.4Mo. Such post- 
interaction binaries are known to exist (90, 95, 96). Owing to its increase in mass, the accretor 
evolves faster, finishing its main sequence and initiating a phase of inverse mass transfer before 
the stripped star depletes core helium. Due to the extreme mass ratio, the mass-transfer is 
unstable, and the system evolves into a common-envelope phase at an age of 130 Myrs. 

The outcome of common-envelope phases is poorly understood in binary evolution (97). 
One-dimensional models can be used to assess whether a merger would occur, or whether the 
hydrogen envelope would be ejected before the helium cores would coalesce. We follow an 
energy prescription for common envelope evolution (44). The binding energy of the envelope 
down to a given mass coordinate m, Eping(m), can be computed as 

Beans f J (- CL 1) dm, (S3) 
m \ rm) 


where M is the total mass of the star, r(m) is the radius of a shell enclosing the mass m, G 


is the gravitational constant, and u is the specific internal energy of the star, for which we also 
consider the energy from hydrogen and helium recombination. The source of energy to unbind 
this material is the change in orbital energy. To estimate the orbital energy released through a 
merger we consider that all layers above m are ejected from the donor, and compute the orbital 


separation anror at which the stripped star (with a radius of 0.18 at the onset of common- 
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envelope) would overfill its Roche lobe. The energy released through the merger, Ego (m), 


is computed then from the orbital energy at this separation, 


Gm Mgtripped 0.18 R 
Eier emj = — ——ÉÓ a m) = © (S4) 
l ) 2arLor(m) RLOP( ) F(Mstrippea/ M) 
where Mstrippea is the mass of the original primary after it has been stripped and 
0.494? 
fla) (S5) 


— 0.63 + In (1 + qi/3) 
is an approximation of the Roche radius (98). The layers that are ejected through the common 
envelope phase are then determined from the condition Epina(Mm) = AcEEmerger(™m), where 
ace € |0, 1] is the efficiency of converting the orbital energy into ejection of the material. 

The binding energy of the donor star and the energy that would be released through a merger 
after ejecting all material above a given mass coordinate are shown in Fig. S11. If orbital energy 
converts with full efficiency into removing the envelope (acg = 1), only 0.08M. of hydrogen 
rich material would remain, indicating that the system would eject the common envelope, halt- 
ing the in-spiral and merging process (Fig. S11). However, the actual efficiency is expected 
to be lower (99). Taking ace = 0.5 or 0.25 the system would retain 0.3Mo and 0.62M, of 
hydrogen rich material, respectively. We perform a simple model of the post-merger evolution 
by constructing a star with a composition profile that combines that of the stripped star with the 
layers that are not ejected from the donor star (Fig. S11). The post-merger evolution depicted 
in Fig. 5 shows the evolution after the merged object reaches a state of thermal equilibrium as a 
core-helium burning object, until it reaches core helium depletion. 

For the merger model with acg = 0.25, we find the best match to the properties of the Wolf- 
Rayet component of HD 45166 occurs at has an age (starting from the initial binary model) of 
133 Myrs, consistent (within 1c) with the age constraint provided by the B7 V component. The 
surface composition is rich in hydrogen, but due to the CNO cycle it has a very low carbon abun- 


dance (Xo = 6.5 x 1075), which is inconsistent with previous measurements of the system (19). 
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That analysis used a wind model that does not describe the magnetosphere (see above). How- 
ever, the observed strength of carbon emission lines (Figs. 1 and 2) in the spectrum indicates that 
the atmosphere of the Wolf-Rayet component contains a significant amount of carbon, which is 
inconsistent with our merger model. We argue that fallback from some of the outermost layers 
of the donor star could have contaminated the merger product with material that has not been 
CNO processed, leading to higher oxygen and carbon abundances observed at the surface in 
the present day. While modeling the merger product beyond core-helium depletion is subject 
to various uncertainties [such as mass-loss in the presence of the magnetic field, (/00)], in our 
model, we find it forms a 1.29 M; O-Ne-Mg core with a 1.04 M. envelope composed primarily 
of helium and with a small amount (0.04M5) of hydrogen. We do not model the evolution of 
the system beyond this point, but we expect the core would continue to grow from shell burning. 
We predict that the Wolf-Rayet component will become a supergiant (R ~ 300 Ro, well within 
its Roche lobe) and explode as a type Ib or IIb electron-capture supernova [ (701), their table 1], 
depending on the amount of hydrogen retained (702, 103). The neutron star remnant is expected 


to possess a magnetic field strong enough to appear as a magnetar (see above). 


Supplementary Text 


Parallels have been drawn between the spectrum of HD 45166 and Wolf-Rayet stars belonging 
to the class of V Sagittae stars (704). This class comprises a few stars that portray Wolf- 
Rayet spectra with reported photometric and spectroscopic variability periods of the order of 
hours that are attributed to orbital periods. Specifically, WR 7a and WR 46 were proposed to 
exhibit properties similar to HD 45166 (/05, 106). However, both WR 7a and WR 46 exhibit 
broad, strong emission lines that resemble those of classical Wolf-Rayet stars, Moreover, no firm 
detection of a magnetic field could be obtained for WR 46 (107). Hence, while an investigation 


of the magnetic properties of WR 7a might still be worth pursuing, neither WR 7a nor WR 46 
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appear to be suitable candidates for massive magnetic helium stars. 

An intriguing similarity is observed when comparing the spectrum of the Wolf-Rayet com- 
ponent to a subset of central stars of planetary nebulae, classified as [WN] or [WC] stars (708). 
The spectral appearance and luminosities of some [WC] and [WN] stars matches well with those 
of the Wolf-Rayet component in HD 45166. However, a distinct difference is that HD 45166 
has no nebula associated with it (58): no evidence for nebulosity is seen in the spectrum, nor is 
it seen in infrared images. While central stars of planetary nebulae are thought to represent the 
final evolutionary stages of solar-type stars following a "standard" post-asymptotic giant branch 
evolution, it is possible that a few impostors are hidden in this population. In the context of the 
evolutionary scenario proposed here, their nebulae could be the remnant of a violent merger 
event, implying that they are "freshly born" magnetic helium stars. Specifically, a few such 
objects were identified that share a similar spectral morphology to HD 45166, and notably the 
Of?p emission complex in the range 4630-4660 A. Examples include PB 8 (108) and Hen 2-108. 
However, a campaign to measure the magnetic fields of a sample of central stars of planetary 
nebulae, including Hen 2-108, did not find evidence for magnetism in these objects (709). It 


remains to be seen whether other [WN] or [WC] stars exhibit strong magnetic fields. 
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Supplementary Figures 


Normalized fluxes 


-50 o 50 100 
Doppler shift [km/s] 


Figure S1: LSD profile of the B7 V companion. The lower (black) curve is the mean Stokes 7 
profile. The upper (red) curve is the measured net circular polarization (Stokes V). The middle 
(blue) curve is the diagnostic null N. No signal is detected inside or outside the line range in 
either V or N. The vertical dashed lines indicate the integration range used for measurement of 
the detection probability and longitudinal magnetic field. 
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Figure 52: Comparison between HD 45166 and NGC 1624-2. Shown are Stokes 
I and V line profiles (see labels in panels) in spectra of HD 45166 (red/solid) and 
NGC 1624-2 (black/dashed) for He 11 \ 4686 (A), Ho +He 11 46560 (B), He 1A 5876 (C), and 
C IV AA 5801,5811 (D). 
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Figure S3: Modeling of the spectropolarimetric signature. Shown are Stokes / and V pro- 
files of O v 45114 (black) compared to LTE polarized spectrum synthesis calculations (red) 
computed assuming a pole-on magnetic dipole of polar strength 45 kG. Such a field geometry 
roughly reproduces the observed profiles. 
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Figure S4: Spectral analysis of the B7 V component. Shown is a comparison between our 
composite spectral model of HD 45166 (red) and the co-added ESPaDOnS spectrum (blue), 
focusing on absorption lines belonging to the B-type component. The composite model in- 
cludes both the Wolf-Rayet and B7 V components, weighed with their respective wavelength- 
dependent light ratios. The spectral lines shown are formed only in the B7 V component, with 
the possible exception of He I 45048 (D), which may be contaminated by the magnetosphere of 
the Wolf-Rayet component. Horizontal black lines are continuum levels. 
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Figure S5: Spectral analysis of the Wolf-Rayet component. Shown is a comparison between 
the co-added IUE spectrum (blue), normalized using the model continua, and three composite 
POWR models calculated with T, = 45,56, and 70 KK (in black, red, and green, respectively) 
for the Wolf-Rayet component. The panels show regions dominated by Fe VI, V, and IV (panels 
A, B, and C, respectively). The contribution of the B7 V component is negligible (see Fig. S6). 
The effective temperature of the Wolf-Rayet component is estimated at T} = 56.0 + 5.0kK. 
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Figure S6: Spectral energy distribution of HD 45166. The observed SED (blue line for IUE 
and FUSE spectroscopy and blue squares for photometry) is compared to the sum (solid red 
line) of the synthetic SEDs of our models for the Wolf-Rayet component (solid black line) and 
the B-type component (dashed green line). The Wolf-Rayet component fully dominates the UV 
flux. The flux excess observed in the infrared could originate in free-free emission stemming 
from the trapped outflow. 
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Figure 57: Frequency analysis of the TESS light curve of HD 45166. A and B: TESS light 
curve from sectors 6 (A) and 33 (B) of HD 45166, with relative brightness in units of milli- 
magnitude (mmag) in the TESS passband ATp on the ordinate axis and time (in BJD) on the 
abscissa. C and D: Amplitude spectra of sector 6 (C) and 33 (D) TESS light curves, with sig- 
nificant frequencies (i.e. exceeding the significance threshold S/N > 4; see text for detailed) 
marked by vertical gray lines. 
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Figure S8: Overview of RVs measured for the two components of HD 45166. (A) Black and 
green symbols show the 199 and 96 RVs derived for the Wolf-Rayet and B7 V components, 
respectively. Each instrument is marked by a different symbol (see legend). The relative RVs of 
the Wolf-Rayet component were vertically shifted to have the same mean as the absolute RVs 
of the B7 V component. (B and C) Short-term variability of both components. The origin of the 
~ 5 hr period of the Wolf-Rayet is uncertain, but could be related to pulsations. We argue (see 
above) that the ~ 1.6d period is due to pulsations in the B7 V component. There is a long-term 
anti-phase motion of the two components, indicating a period of ~ 20 years. Error bars depict 
formal 1c uncertainties (the uncertainties of HERMES and ESPaDOnS RVs are comparable to 
symbol sizes). 
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Figure S9: Spectroscopic signature of gravity-mode pulsations in the B7 V component. 
Shown are co-added nightly ESPaDOnS spectra taken during 2022 Feb 17 and 19 Feb 2022, 
focusing on the Si r1 44128, 4131 doublet of the B7 V component. There are line profile vari- 
ations, with their cores shifting in apparent RV, but their wings being static. This behavior is 
typical for non-radial gravity-mode pulsations (97). 
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Figure S10: Derivation of the mass ratio. The RVs of the B7 V component (y-axis) are shown 
as a function of those of the Wolf-Rayet component (x-axis). Measurements are the same as 
those used for the orbital analysis. The line is a linear regression to the data; its negative slope 
(mass ratio) and intercept are q = 0.60 + 0.13 and C = 20.3 + 2.3. 
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Figure S11: Abundance profiles and envelope binding energy of the pre- and post-merger 
product. (A) abundance profile of the stripped star immediately before the merger. (B, C) 
Abundance profile (B) and binding energy (C) of the donor star at the onset of the in-spiral, 
prior to envelope ejection. Dashed teal curves in panel C show the energy that would be re- 
leased through a merger after ejecting all material above a given mass coordinate for different 
efficiencies acg. The intersection between these curves with Fyina(m), denoted with dashed 
vertical lines, determines the point to which the donor star would be stripped. (D) abundance 
profile of the merger product under the assumption acg = 0.25. Composition profile is con- 
structed from both stars starting from the center with regions depleted in hydrogen in order of 
increasing helium abundance, and then layers with increasing hydrogen abundance. 
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Supplementary Tables 


Table S1: Derived parameters for HD 45166. Provided are the distance, reddening Ep. v, 
age, orbital elements (period P, time of periastron tọ, eccentricity e, argument of periastron w, 
systemic velocity vo, semi-major axis a, mass ratio q, orbital inclination ?, and RV amplitudes 
Kı, K2), mean magnetic-field modulus and line-of-sight components (5) and (B,), minimum 
masses M sin? i, masses M, relative contributions to the light in the V-band I(V), effective 
surface temperatures T., luminosities log L, radii R,, projected rotational velocities v sin 2, ro- 
tational period Pot, and elemental abundances of He, C, N, and O (Xue, Xc, Xn , Xo) in mass 
fractions. Uncertainties are 6896 confidence interval (10). 


Parameter unit B7V Wolf-Rayet 
Distance (adopted) parsec 99135 (110) 

Ep.vy magnitude 0.210 + 0.010 

age Myr 105 + 35 

P day 8200 + 190 

to (MJD) - 49820 + 360 

e - 0.46 + 0.18 

wW degree 132+ 11 

Vo kms"! 5.7 2 1.8 

a au 10.5 + 1.8 

q = Mqwn/ Ms - 0.60 + 0.13 

i degree 49+ 11 

K kms 5.8 t 1.3 9.9 t 1.6 
(B) kG - 43.0 + 2.5 
(Bz) kG < 0.02 13.5 + 2.5 
M sin? i Mo 1.50 +0.74 0.85+0.45 
M Mo 3.40 0.06 2.03 + 0.44 
I(V) - 0.460 + 0.050 0.560 + 0.050 
T; kK 13.00 0.50 56.0 + 6.0 
log L/ Lo - 2.250 + 0.050 3.830 + 0.050 
R, Ro 2.63 + 0.41 0.88 + 0.16 
vsint kms! <10 <10 
Pot day - 124.82 + 0.21 
X ye (adopted) - 0.25 (76) 0.67 (19) 
Xc (adopted) - 0.0024 (76) 0.0059 (19) 
XN (adopted) - 0.00069 (76) 0.0020 (19) 
Xo (adopted) - 0.0057 (76) 0.0015 (19) 
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Table S2: Significant frequencies extracted from the TESS light curves of HD 45166. Fre- 
quencies exceeding the significance criterion (25) in sectors 6 (upper part) and 33 (lower part) 
of the TESS light curve are listed in the first column. The amplitudes and phases are given in 
the second and third columns (Eq. S1). 


Frequency Amplitude Phase 
(d-') (mmag) (radian) 
Sector 6: 
0.6006 + 0.0012 4.23 +0.20 1.618 + 0.094 
0.6825 + 0.0014 3.75 +0.20 —0.11 +0.11 
Sector 33: 


0.59401 + 0.00056 4.45 +0.11 1.114 Æ 0.045 
0.67593 + 0.00057 4.46 +0.12 1.622 + 0.045 
0.75981 + 0.00054 4.68 +0.12 0.142 € 0.043 


Figure S12: * 
Data S1: Derived RVs for the two components of HD 45166. The 1% column lists the instru- 
ment, the 2?4 column lists the MJD, and 3'4, 4^, 5t”, ang 6*^ columns list the RV and the lo 
uncertainty (in kms!) for the Wolf-Rayet and B7 V components, respectively. 
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